In [1]:
import sys
# path to the src folder of the repository
sys.path.append('../src')
import trimesh
import numpy as np
from organ import Organ
from tissue import TissueBlock
from pipeline import Pipeline
from steps import normalize_rigid
from utils.conversions import to_array
from copy import deepcopy
from tqdm.auto import tqdm
Load a mesh as a point cloud¶
In [2]:
# load source and reference organs
# source_atlas = Organ(path='/Users/dhruviljoshi/Documents/koding/cns/hra-amap/data/Ovary/Source/generic_ovary_updated.glb')
# hra_atlas = Organ(path='/Users/dhruviljoshi/Documents/koding/cns/hra-amap/data/Ovary/Reference/VH_F_Ovary_L.glb')
source_atlas = Organ(path='../data/3d-vh-m-kidney-l.glb')
hra_atlas = Organ(path='../data/3d-vh-m-kidney-r.glb')
Visualize point cloud (before)¶
In [3]:
# normalize (optional)
output = normalize_rigid(source=deepcopy(source_atlas.pointcloud), target=deepcopy(hra_atlas.pointcloud))
normalized_source_atlas = to_array(output.output['Source'])
normalized_target_atlas = to_array(output.output['Target'])
In [4]:
source_pc = trimesh.PointCloud(normalized_source_atlas, colors=np.tile(np.array([255, 0, 0, 1]), (len(normalized_source_atlas), 1)))
hra_pc = trimesh.PointCloud(normalized_target_atlas, colors=np.tile(np.array([0, 0, 255, 1]), (len(normalized_target_atlas), 1)))
before_scene = trimesh.Scene([source_pc, hra_pc])
# before_scene.export('/Users/bhargavdesai/Desktop/test/before.glb')
before_scene.show()
Out[4]:
In [5]:
trimesh.Scene([source_pc]).show()
Out[5]:
In [6]:
trimesh.Scene([hra_pc]).show()
Out[6]:
Register point cloud¶
In [7]:
# instantiate the registration pipeline
pipeline = Pipeline(name='Base Registration', description='Base Registration', params='../configs/params.yaml')
In [8]:
# run registration
projections = pipeline.run(source=source_atlas, target=hra_atlas)
Visualize point cloud (after)¶
In [9]:
projected_pc = trimesh.PointCloud(projections.registration.vertices, colors=np.tile(np.array([255, 0, 0, 1]), (len(projections.registration.vertices), 1)))
hra_pc = trimesh.PointCloud(hra_atlas.vertices, colors=np.tile(np.array([0, 0, 255, 1]), (len(hra_atlas.vertices), 1)))
after_scene = trimesh.Scene([projected_pc, hra_pc])
# after_scene.export('/Users/bhargavdesai/Desktop/AMap/Figures/Figure 1/after.glb')
after_scene.show()
Out[9]:
In [10]:
projections.registration.show()
Out[10]:
Heatmap¶
In [11]:
from sklearn.preprocessing import minmax_scale
In [12]:
sampled = projections.registration.simplify_quadric_decimation(0.5)
In [13]:
print(hra_atlas.is_watertight)
print(projections.registration.is_watertight)
False False
In [14]:
print(projections.registration.vertices)
[[-0.06591257 0.23947445 0.01211492] [-0.06441004 0.23920059 0.0127866 ] [-0.06627862 0.24042541 0.01224176] ... [-0.07925211 0.24227494 0.00386304] [-0.0787062 0.24148023 0.00421847] [-0.07894563 0.24225689 0.00393145]]
In [15]:
print(type(hra_atlas))
print(type(projections.registration.vertices))
<class 'organ.Organ'> <class 'trimesh.caching.TrackedArray'>
In [16]:
# sd = trimesh.proximity.signed_distance(hra_atlas, projections.registration.vertices)
In [17]:
sd = trimesh.proximity.signed_distance(hra_atlas, sampled.vertices)
In [18]:
# rescale
sd_scaled = minmax_scale(sd, feature_range=(-1, 1), axis=0, copy=True)
# create colors
cmap = 'jet'
colors = trimesh.visual.interpolate(sd, color_map=cmap)
# rebuild the registered organ
# heatmap = trimesh.Trimesh(vertices=np.array(projections.registration.vertices), faces=np.array(projections.registration.faces), vertex_colors=colors)
heatmap = trimesh.Trimesh(vertices=np.array(sampled.vertices), faces=np.array(sampled.faces), vertex_colors=colors)
In [19]:
# heatmap.export('/Users/bhargavdesai/Desktop/AMap/Figures/Figure 1/errormap.glb')
heatmap.show()
Out[19]:
Histogram¶
In [20]:
import seaborn as sns
In [21]:
sns.set_style('darkgrid')
sns.set_theme('paper')
In [22]:
hist = sns.histplot(data=sd_scaled).set_xlabel('Error (Signed Distance)')
fig = hist.get_figure()
# fig.savefig("/Users/bhargavdesai/Desktop/AMap/Figures/Figure 1/hist.png")
In [ ]:
!jupyter nbconvert --to html RegistrationErrorVisualization.ipynb
In [ ]: